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While landmark studies iiave sliown tliat mlcrobiota activate and educate liost Immunity, iiow immune systems sliape 
microbiomes and contribute to disease is incompletely characterized. Primary immunodeficiency (PID] patients suffer 
recurrent microbial infections, providing a unique opportunity to address this issue. To investigate the potential influence 
of host immunity on the skin microbiome, we examined skin microbiomes in patients with rare monogenic PIDs: hyper-lgE 
{STAT3-deficlent), Wiskott-Aldrich, and dedicator of cytokinesis 8 syndromes. While specific immunologic defects differ, 
a shared hallmark is atopic dermatitis (AD)-like eczema. We compared bacterial and fungal skin microbiomes (41 PID, 13 
AD, 49 healthy controls] at four clinically relevant sites representing the major skin microenvironments. PID skin dis- 
played increased ecological permissiveness with altered population structures, decreased site specificity and temporal 
stability, and colonization with microbial species not observed in controls. Including Clostridium species and Serratia 
marcescens. Elevated fungal diversity and increased representation of opportunistic fungi ^Candida, Aspergillus) supported 
increased PID skin permissiveness, suggesting that skin may serve as a reservoir for the recurrent fungal Infections observed 
in these patients. The overarching theme of increased ecological permissiveness in PID skin was counterbalanced by 
the maintenance of a phylum barrier in which colonization remained restricted to typical human-associated phyla. 
Clinical parameters. Including markers of disease severity, were positively correlated with prevalence of Staphylococcus, 
Corynebacterium, and other less abundant taxa. This study examines differences in microbial colonization and community 
stability in PID skin and Informs our understanding of host-microbiome Interactions, suggesting a bidirectional dialogue 
between skin commensals and the host organism. 



[Supplemental material is available for this article.] 

Although the relationships between mlcrobiota, host, and disease 
are complex. Investigators have long relied upon Koch's postu- 
lates to impute causation between microorganism and disease. 

However, the requirement that disease-causing microorganisms 
uniformly recapitulate disease in healthy individuals was com- 
promised when Koch discovered asymptomatic carriers of Vibrio 
cholera and Salmonella typhi, documenting the Important dis- 
tinction between apparent colonization and infection. Poten- 
tially pathogenic bacteria that colonize humans are influenced 
not only by host factors, but also by their native microbial com- 
munity. This paradigm requires that diseases of microbial origin 
must be investigated within the context of their microbial com- 
munity, host factors, and immunity. 

Genome-wide surveys of the microbes inhabiting the differ- 
ent microenvironments of the human body (e.g., gut, skin, oral) 
(The Human Microbiome Project Consortium 2012), model or- 
ganisms (Chung et al. 2012), and geographic regions (Yatsunenko 
et al. 2012) have suggested the presence of strong selective control 
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on the microbial communities that can colonize and persist in 
particular environments. The identity and relative representation 
of microbes at a given site, the overall microbial density (bio- 
burden) and diversity are all components that can be constrained 
by extrinsic or intrinsic factors, such as nutrient accessibility, 
temperature, oxygen, or products produced by other microbes. 

Innovations in sequencing technology have enabled explo- 
rations of the complexity of the vast human-associated mlcro- 
biota. Recent reports have examined shifts in microbial species and 
communities associated with diseases as varied as atopic dermati- 
tis, inflammatory bowel disease, colorectal cancer, obesity, and 
metabolic syndrome (Baiter 2012). Microbes also promote human 
health by performing vital functions such as aiding in digestion 
and educating the Immune system (Nalk et al. 2012; Upadhyay 
et al. 2012). Less well understood is the role that human immunity 
plays in shaping resident microbial communities. Some environ- 
mental microbial communities are extremely simple, containing 
just a few species (e.g., acid mine drainage biofilms) (Denef et al. 
2010), while others are enormously complex including hundreds 
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of phyla, Including the Guerrero Negro hypersaUne mats (Harris 
et al. 2013), soil, or ocean. Animal-associated microbial symbioses 
range from monospecific association of Vibrio fischeri inhabiting 
the light organ of Euprymna scolopes (squid) to the modestly com- 
plex communities inhabiting the human body (Dethlefsen et al. 
2007). Even at its most complex, the human system appears to be 
colonized by approximately 10 phyla, an order of magnitude fewer 
than many environmental ecosystems. An open question is the 
extent to which symbiotic bacteria are selected for mutualistic 
interactions with the human host based on the limited ecological 
niches of human tissues versus a more active modulation of mi- 
crobial communities by the immune system. 

The microbiome of healthy human skin shows significant 
topographiccil variation between sebaceous, dry, and moist micro- 
environments. Likewise, the temporal stability of bacterial com- 
munity membership and structure varies by skin site (Costello et al. 
2009; Grice et al. 2009; The Human Microbiome Project Consor- 
tium 2012). These attributes differ significantly in disease states 
such as atopic dermatitis. Skin disease severity in atopic dermatitis 
patients was closely associated with marked shifts in the composi- 
tion, diversity, and interpersonal and temporal variation of the 
resident microbial communities (Kong et al. 2012). 

Studying primary immunodeficiency (PID) patients provides 
an opportunity to evaluate the degree to which altered immunity 
may influence the human microbiome and how, in turn, micro- 
biota may interact with the host to cause disease. As cutaneous 
immunity has been shown to be distinct from that of other organs 
Oiang et al. 2012; Naik et al. 2012), we surveyed the skin micro- 
biomes of three different groups of PID patients whose genetically 
defined diseases are characterized by AD-like skin disease, di- 
minished memory T and B cells, and variable eosinophilia and 
elevated IgE. Importantly, despite this broad phenotypic overlap, 
each of these diseases has its own distinct monogenic mutations 
and other extracutaneous clinical manifestations: Signal trans- 
ducer and activator of transcription 3 (STAT3) mutations cause 
autosomal dominant hyper-IgE syndrome (STAT3-HIES; Job syn- 
drome, OMIM 243700); mutations of dedicator of cj^okinesis 8 
(DOCKS) cause autosomal recessive hyper IgE syndrome (DOCKS 
deficiency, OMIM 611432); and mutations in WAS cause X-linked 
recessive Wiskott-Aldrich syndrome (WAS; OMIM 301000). 

STAT3-HIES arises from dominant-negative STAT3 mutations 
that prevent Thl7 cell differentiation, abrogating CD4-I- inter- 
leukin (1L)-17 production, which is important in epithelial fungal 
and bacterial infections as well as defensin expression (Milner et al. 
2010). Clinical manifestations include recurrent staphylococcal 
infections, eczema, sinopulmonary infections, mucocutaneous 
candidiasis, and abnormalities in connective tissue, vessels, teeth, 
and bones, as well as an increased rate of lymphoma (Holland et al. 
2007; Minegishi et al. 2007). DOCKS deficiency more generally 
impairs T-cell differentiation, including Thl7 cells (Engelhardt 
et al. 2009; Zhang et al. 2009), and is postulated to result from 
cytoskeletal abnormalities that may inhibit proper migration and 
differentiation of immune cells. In addition to recurrent sino- 
pulmonary infections and staphylococcal skin infections, DOCKS 
deficiency is associated with eczema, persistent viral infections of 
the skin, and a predisposition to malignancies (Renner et al. 2004). 
WAS results from X-linked recessive mutations in WAS. WAS has 
cytoskeletal abnormalities that lead to defects in T- and B-lym- 
phocyte function via deficiencies in locomotion, signaling, and 
immune synapse formation (Ochs and Thrasher 2006). Clinical 
manifestations of WAS include eczema, thrombocytopenia, neu- 
tropenia, and a predisposition to autoimmunity and malignancy. 



PID patients suffer from recurrent infections caused by readily 
isolated common bacteria and fungi (e.g.. Staphylococcus aureus and 
Candida albicans) as well as uncommon opportunistic pathogens. 
Diagnosis typically relies on culture-based analyses; however, cul- 
tivation methods are limited in their ability to investigate complex 
microbial communities or to detect fastidious microbiota. Culture- 
independent metagenomic approaches, e.g., 16S ribosomal RNA 
(rRNA) sequencing, greatly increase the range of detectable mi- 
crobes and allow assessment of the relative abundance of members 
of the microbial community. The 16S rRNA gene is present in all 
bacteria/archaea and contains both variable regions, which enable 
taxonomic classification, and conserved regions, which serve as 
binding sites for PGR primers. 

Here, we present the results of a culture-independent com- 
parative analysis, using 16S rRNA sequencing, of the skin micro- 
biota of three groups of PID patients compared with classical AD 
patients and healthy controls. Our goal was to examine alterations 
in skin microbial communities in human immunodeficiencies as 
a basis for understanding the forces exerted by both the ecological 
niche and immune selection. 

Results 

Data generated for microbiome analysis 

To investigate the extent to which the skin microbiome is altered 
in PID patients compared with controls, we analyzed a total of 
170,167 full-length 16S rRNA Sanger sequences (mean 60S reads/ 
sample; greater than 199 reads minimum, 2S0 samples), 6,S50,352 
Vl-3 (average 11,650 reads/sample; greater than 1000 reads min- 
imum, 595 samples), and 1,496,770 internal transcribed sequence 
regions 1 (ITSl; 13,607 mean reads/sample; greater than 1051 
reads minimum, 90 samples) pyrosequences for four skin/nares 
sites in 49 healthy controls, 13 classical AD, 25 STAT3-HIES, 10 
WAS, and six DOCKS patients (Supplemental Table SI). Patient 
characteristics are described in Table 1. In general, the correlations 
of relative abundances of 17 major skin taxa obtained via Sanger 
compared with 454 sequencing were robust (partial Spearman 
correlation adjusting for multiple measurements 0.42 < p < 0.81; 
2.6 X 10 ' < P < 9.9 X 10"^^^). Our findings justified use of se- 
quences from the two platforms interchangeably for taxonomy- 
based analysis, while noting that full-length Sanger sequencing 
differentially recovered sequences for Enterobacteriaceae (Supple- 
mental Fig. SI; Supplemental Table S2). For analyses based on 
operational taxonomic units (OTU; 97% similarity cutoff), such as 
alpha and beta diversity metrics, only 454 data were used. Rare- 
faction curves indicated that sampling for both Sanger and 454 
sequencing provided sufficient coverage to analyze the dominant 
members of the bacterial communities (Supplemental Fig. S2; 
Supplemental Table S3). 

Since PID patients often use systemic and/or topical medica- 
tions to prevent and/or treat infections, we included both healthy 
volunteers and AD patients undergoing treatment as controls 
(Supplemental Table S4). For each patient cohort, we examined the 
use of antimicrobial agents and immunosuppressants for possible 
associations with changes in the skin and mucosal bacterial com- 
munities. Bacterial community diversity was measured by the 
Shannon diversity index, an ecological measure of microbial 
communities that considers the following: (1) richness, or the total 
number of bacterial types, and (2) evenness, or the relative pro- 
portion of these bacterial types. Within each PID group, patients 
who had received systemic treatment within 30 d (antimicrobials. 
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Table 1. Summary of patient characteristics 



Controls Primary immunodeficiency 



Croup 




Healthy 


Atopic dermatitis 


Hyper-lgE 


Wiskott-Aldric 


DOCKS deficiency 


Subjects analyzed 




49 


13 


25 


10 


6 


Age, median (range) 




20 (2-40) 


7 (2-1 7) 


1 3 (4-37) 


1 6 (2-27) 


21 (8-26) 


Female:male 




19:30 


5:8 


11:14 


0:10 


3:3 


Sites analyzed 




Af, Pf, Vf, Ra, N 


Af, Pf, Vf, N 


Af, Vf, Ra, N 


Af, Vf, Ra, N 


Af, Vf, Ra, Pf, N 


Objective SCORAD/ 


median (range) 


0 


28 (7.9-60.1) 


22 (0-69) 


6 (0-18.6) 


11.1 (0-25) 


IgE in lU/mL, median 


(range) 


b 


51.4 (12.9-2790) 


11510 (242-67,122) 


37.6 (12.3-665) 


8765 (279-37,822) 


Absolute eosinophils in K/p.L, median (range) 


c 


0.64 (0.34-0.69) 


0.735 (0.03-2.20) 


0.32 (0.03-0.64) 


0.215 (0.08-0.33) 


Absolute lymphocyte 


in K/^JLL, median (range) 


d 


2.1 (2.10-3.70) 


2.73 (1 .42-4.93) 


2.74 (1.25-4.72) 


0.645 (0.41-1.31) 



^Scoring atopic dermatitis. 

"Normal range IgE: 0-90 lU/mL. 

^Normal range absolute eosinophils (K/^JLL): 0.04-0.54. 

''Normal range absolute lymphocytes (K/|xL): 1.18-3.74. 



steroids) and/or topical treatment within 7 d (antimicrobials, ste- 
roids) did not harbor generalized changes in bacterial community 
diversity when compared to PID patients not receiving recent 
treatment (Supplemental Fig. S3; Supplemental Table S4). We did 
not observe an obvious shift in the bacterial community diversity 
related to recent treatment, distinct from findings identified in AD 
patients in which use of AD treatments was associated with an 
increase in bacterial community diversity compared with an un- 
treated state of AD disease exacerbation (Kong et al. 2012). In- 
terestingly, a recent study demonstrated that the skin microbiome 
in mice did not change after oral antibiotic treatment, in contrast 
to the gut microbiome (Naik et al. 2012). Finally, as skin microbiota 
differ based on age (Oh et al. 2012), patients and healthy controls 
were matched based on the two major groupings of prepubertal/ 
early-pubertal versus late adolescent/adult (Supplemental Table 
S6), with the exception of the retroauricular crease (Ra), for which 
no prepubertal/early-pubertal data existed. 

Given the site specificity observed in the skin microbiome 
(Costello et al. 2009; Grice et al. 2009; The Human Microbiome 
Project Consortium 2012) and the tendency for dermatoses to 
present in characteristic sites of predilection, we examined multi- 
ple skin sites in our patients: the ante- 
cubital fossae (Af) and the Ra as sites of 
disease predilection; volar forearm (Vf) as 
an adjacent control skin site that is not 
characteristically affected; and the nares 
(N), a potential reservoir for pathogens 
such as S. aureus and Streptococcus pneu- 
moniae (von Eiff et al. 2001). For AD, 
STAT3-HIES, WAS, and DOCKS patients, 
skin disease severity was assessed quanti- 
tatively with objective SCORAD (scoring 
AD), a well-validated clinical assessment 
tool (Williams et al. 1994, 1996; Kunz 
et al. 1997; Oranje et al. 2007). Laboratory 
studies in PID patients included serum IgE 
and a complete blood count (Supplemen- 
tal Table S4). Representative images com- 
paring skin disease manifestations in each 
patient group are shown in Figure 1, A 
through D. Because our largest PID patient 
cohort was the STAT3-HIES group, with 



between STAT3-H1ES patients and healthy controls, with addi- 
tional analyses including the other PID/AD cohorts to broaden 
insights into the microbiomes of these disorders. 

Altered permissiveness to microbial colonization in PID 
patients 

We sought to determine whether underlying immunodeficiency 
altered the bacterial communities present on human skin, specif- 
ically to include non-native taxa or unique community structures 
on PID skin. We classified the sequence data to (1) determine if 
unique bacteria colonized the skin of PID patients and (2) assess 
the taxonomic differences between PID patients and healthy 
controls (Fig. 2A,B; Supplemental Fig. S4A-D; Supplemental Table 
S9). We analyzed sites of disease predilection (Af, Ra) and observed 
that the environmental microbe and opportunistic pathogen 
Proteobacteria Serratia marcescens was increased at both sites in 
STAT3-HIES patients, as well as the adjacent arm site (Af, Vf: ad- 
justed P < 0.05; Ra: N.S.) (Fig. 2A). Although S. marcescens is 
ubiquitous in the environment and is commonly found in 
households, this bacterium has not previously been identified as 
a component of the skin microbiome of healthy individuals or AD 




ljUj 



numerous recurring visits, we focused 
the bulk of our analyses on comparisons 



Figure 1. Representative clinical images of disease severity in the different patient groups. {A) Non- 
PID atopic dermatitis, (fi) Hyper IgE syndrome, (C) Wiskott-Aldrich, and (D) DOCKS deficiency. 
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Figure 2. Bacterial taxonomic classifications show colonization with unique taxa and altered repre- 
sentation of diverse taxa. Relative abundances of 14 major phyla-family taxonomies in the antecubital 
fossa (Af) are shown in A and retroauricular crease (Ra) in 6. Shown are seven representative healthy 
controls and all Sanger-sequenced STAT3-HIES patients. Full versions of all skin and nares sites are in 
Supplemental Figure S4. (C) Inlaid pie charts compare the mean relative abundances of major skin 
genera across patient groups. Mean ± SEM for C are shown in Supplemental Figure S5 and Supple- 
mental Tables S9 and SI 0. Patient identifiers are in Supplemental Table SI 3. 



patients, independent of treatment with antimicrobials. S. mar- 
cescens is a nosocomial pathogen associated with wounds, urinary 
tract infections, and bacteremia, but it is not known to cause 
clinical disease in STAT3-H1ES patients. Moreover, S. marcescens 
was observed infrequently in the other two PID groups (Supple- 
mental Fig. S4) but can cause disease in some patient populations, 
e.g., chronic granulomatous disease. 

The Af (Fig. 2A) and Vf (Supplemental Fig. S5) are typified by 
intermediate bacterial diversity with moderate amounts of Pro- 
pionibacterium, Staphylococcus, and Corynebacterium, while the se- 
baceous Ra (Fig. 2B), which presents a different cutaneous envi- 
ronment, is characterized in healthy individuals by low bacteria 
diversity with a predominance of Propionibacterium. Further anal- 
yses of the skin sites (Table 2; Supplemental Table S9) showed 
a significant overrepresentation of Firmicutes in the STAT3-HIES 
skin, particularly of the Clostridiales class {Anaerococcus, Finegoldia, 
and/or Peptoniphilus, FDR-adjusted P £ 0.1). The common skin 
bacteria Corynebacterium and Staphylococcus were highly over- 
represented on STAT3-H1ES skin compared with healthy controls 
(adjusted P < 0.05) (Fig. 2C). Bacteroidetes (particularly Porphy- 
romonas and/or Cloacibacterium, adjusted P < 0.1 Af, Vf) and the 
abundant skin commensal Propionibacterium (adjusted P = 0.1) 
were significantly depleted at sites of disease predilection (Af, Ra) 



in the STAT3-HIES individuals. With the 
finding of a distinct species {Serratia mar- 
cescens) colonizing STAT3-H1ES skin, 
a possible hypothesis is that a deficiency 
in the host immune system may in- 
crease the skin's ecological permissive- 
ness, e.g., allowing atypical microbiota 
to colonize the skin. However, we did 
not identify a generalized increase in 
the number of different taxa on PID 
skin, even in patients who had not re- 
ceived antimicrobials. Although we 
observed a unique colonization by 
S. marcescens, the bacterial phylum that 
includes S. marcescens, Proteobacteria, 
also includes other species common to 
the human skin. Thus, the ability of 
human skin to host certain bacterial 
phyla remained restricted to those pre- 
viously identified. 

To assess which taxa contributed to 
the variation observed between samples, 
we performed principal coordinates anal- 
ysis (PCoA) using the Yue-Clayton theta 
(6) similarity coefficient (Clayton and Yue 
2005), which takes into consideration the 
number of bacterial species present and 
their relative abundances in the two 
communities being compared. Spearman 
correlation of the relative abundances 
of different taxa to the axes, shown by 
biplot lines, identifies taxa that most 
significantly contribute to axis variation. 
At sites of disease predilection, the Af 
and Ra (Fig. 3; for additional skin and N 
sites, see Supplemental Fig. S6), the pre- 
dominance of Staphylococcus, Propioni- 
bacterium, Corynebacterium, and members 
of the Clostridiales family appeared to be 
the primary drivers of variation between individuals and groups. 
In particular. Staphylococcus, and to a lesser extent, Co- 
rynebacterium, contributed most heavily to the axis variation as- 
sociated with PID and AD patients, while Propionibacterium was 
heavily associated with healthy controls. To quantify the level of 
similarity of bacterial community structures in our patient cohort 
and healthy controls, we calculated and compared interpersonal 
variation: (1) between individuals of a patient group ("within- 
group"), (2) between different patient groups and controls ("vs. 
controls"), and (3) between different PID groups ("between PID") 
(Supplemental Fig. S7). Interpersonal variation between the micro- 
bial communities of the STAT3-H1ES individuals was significantly 
lower than in healthy controls and other PID groups at all sites 
(multiple comparison test after Kmskal-Wallis, P < 0.05) (Supple- 
mental Table S8) except the Ra, which displayed low interpersonal 
variation irrespective of patient group. Interpersonal variation 
within STAT3-H1ES, DOCKS, and WAS groups was generally similar 
(P < 0.05), again with the exception of the Ra. At sites of eczematous 
disease predilection (Af and Ra), STAT3-H1ES and DOCKS in- 
dividuals were significantly more similar to each other than to WAS 
individuals (P < 0.05, Af), while WAS individuals were generally 
significantly more similar to controls (P < 0.05 for Af, Ra, Vf) than 
other PID individuals. 
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Table 2, Select abundant taxonomies differentially abundant between controls and STAT3-HIES individuals 
Classification 



Site Pliylum/Order/Genus/Species CTRL (mean ± SE) STAT3-HIES (mean ± SE) Adjusted P-value 



Actinobacteria 


Actinomycetoles 


Corynebacterium 


9.1 ± 1 .3% 


1 7.2 + 2.6% 


0. 04678 






Propionibacterium 


24.2 ± 3.0% 


10.9 + 2.7% 


0.1 0688 


Bocteroiaetes 


Flovobacteriaies 


Porphyromonas 


r\ Q J. n 1 0/ 
U.y ± U. 1 70 


U.3 ± U. 1 % 


A A/1 Q T 7 






Cloacibacterium 


3.0 ± 0.7% 


0.0 ± 0.0% 


0.04678 


Firmicutes 


Bacilliales 


Staphylococcus 


11.3 ± 1.7% 


27.7 ±4.0% 


0.03092 






S. aureus 


0.4 ± 0.2% 


10.3 ± 3.4% 


0.02281 






S. epidermtdis 


C /I J. 1 AO/ 


1 Q T _i_ T "70/ 
1 D.2 ± 2./% 


A A/1 QO 






S. hoemolyticus 


U.z ± U.Uvo 


1 /I J- n A OA. 

\ A ± U.4vo 


A AAACA 
U.UUUjU 




Clostriaiales 


Anaerococcus 


0.6 ± 0.1 % 


T o 1 r» yi n/ 
1 .O + U.4% 


0.04480 






Finegoldia 


0.2 ± 0.0% 


1 .6 ± 0.4% 


0.00001 






Peptoniphilus 


U.4 ± U. 1 tX) 


U.b ± U.2% 


U.UoD/o 


Proteobactena 


Betaproteobacteria 


Diaphrobacter 


2.1 ± 0.6% 


A A _L A AO/ 
0.0 ± 0.0% 


0.04678 




Cammaproteobacteria 


Serratia 


0.0 ± 0.0% 


8.1 ± 3.3% 


0.00218 


Firmicutes 


Bacilliaies 


Staphylococcus 


lU./ ± 1 .3% 


52.5 ± 4.0% 


A AAAAT 
U.UU002 






S. aureus 


3.2 ± 0.9% 


5.0 + 1 .9% 


0.74406 






S. epidermtdis 


C f\ A- f\ 00/ 

o.u ± u.y% 


24. o ± 4.2% 


A AAAAA 
U.UUUUU 






S. haemolyticus 


0.0 ± 0.0% 


0.5 ± 0.2% 


0.0241 2 






Streptococcus 


6.2 ± 1 .1 % 


5.3 ± 1 .4% 


0.98904 


Actinobacteria 


Actinomycetales 


Corynebacterium 


4.i5 ± 1 .z% 


1 O A _i_ /I CO/ 

1 y.u ± 4.3% 


A AA1 C 

U.UU 1 6i 






Propionibacterium 


67.5 ± 4.8% 


20.8 ± 6.8% 


0.00004 


Firmicutes 


Bacilliales 


Staphylococcus 


1 o.y ± 3.6% 


U "7 Q J. C TO/ 

3/. 8 ± 5.2% 


0.00556 






S. aureus 


U.O ± O.U% 


y. 1 ± 3.9% 


U.UUc524 






S. epidermidis 


1 3.0 ± 2.9% 


24.6 + 5.4% 


0.0421 9 






S. haemolyticus 


0.0 ± 0.0% 


2.1 ± 1 .6% 


0.07403 






Finegoldia 


0.3 ±0.1% 


2.0 ± 0.6% 


0.00197 


Actinobacteria 


Actinomycetales 


Corynebacterium 


7.9 ± 1 .4% 


1 7.6 ± 2.6% 


0.01084 


Bacteroidetes 


Flovobacteriaies 


Cloacibacterium 


2.1 ± 0.6% 


0.0 ± 0.0% 


0.10420 






Porphyromonas 


1.2 ±0.2% 


0.5 ±0.1% 


0.07706 


Firmicutes 


Bacilliales 


Staphylococcus 


6.3 ±1.1% 


18.3 ±2.7% 


0.02686 






S. aureus 


0.6 + 0.3% 


7.5 ± 2.4% 


0.05623 






S. epidermidis 


1 .6 ± 0.2% 


7.1 ±1.5% 


0.00002 






S. haemolyticus 


0.2 ± 0.1% 


1 .4 ± 0.4% 


0.00042 




Clostridiales 


Anaerococcus 


0.5 ± 0.1% 


2.0 ± 0.4% 


0.00001 






Finegoldia 


0.2 ± 0.0% 


2.2 ± 0.5% 


0.00000 






Peptoniphilus 


0.3 ± 0.1% 


1.2 ± 0.3% 


0.00084 


Proteobacteria 


Cammaproteobacteria 


Serratia 


0.0 ± 0.0% 


8.4 ± 3.3% 


0.00473 






Haemophilus 


0.8 ±0.1% 


0.1 ± 0.0% 


0.01084 



In general, consistent with tlie beta diversity metrics siiowing 
more similarity between the WAS and control groups, there were 
few genera whose mean relative abundances were significantly 
different between the WAS patients and controls, irrespective of 
site. While the DOCKS patient population was too small to dem- 
onstrate statistically significant differences between groups, sev- 
eral of the trends observed in the STAT3-HIES population, such as 
elevated Clostridiales and Staphylococcus epidermidis levels, were 
also observed in this population (Supplemental Fig. S6). 

Finally, to identify whether there were key Indicator taxa that 
could be used to classify an individual as belonging to the STATS- 
HIES immunodeficient patient class, we used Random Forests, 
a supervised machine learning technique (Knights et al. 2011). To 
determine those most differentiating taxa for each site, we ranked 
taxa by its mean decrease in accuracy, a measure of its differenti- 
ating power calculated by estimating the increase in classification 
error caused by removing that taxa from the set of predictors. The 
error rate for each site ranged from 1 0% (Vf) to 1 8% (Ra; Af : 1 6%; N: 
13%). Consistent with our over- and underrepresentation analy- 
sis, taxa of the Clostridiales order, P. acnes, S. aureus, Serratia, and 
Cloacibacterium were key indicator taxa differentiating the two 
groups consistently across skin sites, with key differentiating taxa in 



the N S. epidermidis, Enterococcus, and unclassified Actinomycetales 
(Supplemental Table S14). 

Ecological differences in microbial communities of PID patients 

To further explore microbial community changes, we compared 
bacterial community diversities in the skin of the PID cohort and 
controls using Shannon diversity, an ecological measure of com- 
munity richness and evenness. At the Af, diversity in STAT3-HIES 
individuals trended toward being lower compared with healthy 
controls (P = 0.08) (Supplemental Table S7), with levels more 
comparable to those of AD patients (P > 0.05) (Fig. 4A). At char- 
acteristically unaffected sites such as the Vf and N, diversity was 
comparable between all PID patient groups, AD, and healthy 
controls (Kruskal-Wallis test, P > 0.05; STAT3-HIES N excepted P = 
0.003). While the Af, Vf, and N communities are moderately 
complex, the sebaceous Ra is characterized by markedly low bac- 
terial diversity with a predominance of Propionibacterium. Notably, 
we observed a loss of this homogeneity with the bacterial diversity 
at the Ra being significantly elevated in all PID cohorts (Kruskal- 
Wallis test, P - 0.01). This ecological measure of community di- 
versity mirrored results obtained by taxonomic classification of the 
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Figure 3. Taxonomies associated with variation within and between individuals 
ordinates analysis (PCoA) of the Yue-Clayton theta coefficient, which calculates the similarity between 
two samples based on (1) number of species in common between two samples, and (2) their relative 
abundances. Samples that have similar principal coordinates appear closer together, i.e., are more 
similar. Biplot lines indicate the most significant unique consensus taxonomies contributing to variation 
along axis 1 ; Spearman correlations (p) are with axis 1 . Length of biplot lines reflects the contribution of 
that taxa to the top three axes. Associated P-values < 2.2 x 10"^^. Percentage variation attributed to an 
axis is indicated. Antecubital fossa (Af) and retroauricular crease (Ra) are shown in A and 6, respectively; 
nares and volar forearm are shown in Supplemental Figure S6. 



in healthy individuals are often strikingly 
dissimilar from site to site (Costello et al. 
2009; Grice et al. 2009). However, we 
observed significant blunting of this site- 
to-site variation, or decreased site-speci- 
ficity, in the bacterial community struc- 
ture of the DOCKS and STAT3-H1ES 
microbiomes (P = 0.02 and P = 0.0001, 
respectively) and a similar trend in WAS 
(P = 0.07) (Fig. 4B; Supplemental Table 
S8). 

Another trend that might reflect in- 
creased permissiveness of PID skin would 
be aberrant stability of the microbial 
communities over time. In healthy con- 
trols, longitudinal stability of microbial 
communities is highly site-specific, asso- 
ciated most closely with microenviron- 
ment characteristics, including moist or 
sebaceous skin. Since immunodeficiency 
may alter temporal stability of skin mi- 
crobial communities, we compared the 
longitudinal stability between healthy 
controls and STAT3-HIES patients, cal- 
culating the mean theta index across two 
to four timepoints for STAT3-HIES pa- 
tients (Fig. 4C; Supplemental Table S8). 
Across multiple sites, microbiomes of 
STAT3-H1ES patients were significantly 
less stable over time (Af: P = 0.19; N: P = 
0.02; Ra: P = 0.002; Vf: P = 0.02) com- 
pared with healthy controls. Alterations 
in skin bacterial diversity, blunting of 
site-to-site variability, and increase in 
temporal instability, taken together with 
the ability of unique taxa to associate with 
PID skin, suggest that immunodeficiency 
may result in dysregulation of the typical 
host-specific constraints on community 
composition. 



PID skin in which expansion in the number of phyla able to col- 
onize the PID skin was not observed. However, dysbiosis, in terms 
of bacterial community diversity, existed in a site-specific manner, 
i.e., sites of disease predilection. 

While overall community diversity was not significantly in- 
creased in PID skin, we hypothesized that PID could potentially 
increase permissiveness of microbial colonization in the skin, such 
that community composition at different sites is dysregulated. To 
compare the similarities of bacterial community structures be- 
tween different sites, we calculated the mean theta index across 
three diverse microenvironments of the skin (Af, Ra, and Vf : moist, 
sebaceous, and dry, respectively) and the N. The skin microbiomes 



Association of tiie microbiome 
witii clinical disease manifestations 
and severity 

Since skin disease in PID patients is often 
associated with S. aureus, we classified 
Staphylococcus to the species level (Fig. 5A; 
Supplemental Fig. S4; Supplemental Ta- 
ble SIO). While staphylococci were gen- 
erally more abundant in the AD group than in the STAT3-H1ES Af 
and Vf (adjusted P = 0.03, 0.06, respectively), we observed signifi- 
cant differences in the species composition in the STAT3-H1ES and 
DOCKS groups. In contrast to AD patients, whose staphylococci are 
primarily S. aureus (Fig. 5A; Supplemental Figs. S4, S5), at all sites 
S. epidermidis was significantly enriched in the STAT3-H1ES group 
(adjusted P < 0.05), and Staphylococais haertiolyticus, an opportu- 
nistic pathogen, was less abundant but also significantly enriched 
(adjusted P s 0.07). While S. aureus was also significantly enriched 
at skin sites (adjusted P £ 0.06; equally abundant in the N), sur- 
prisingly, it comprised a much smaller proportion of staphylococci 
in the STAT3-H1ES and DOCKS groups as compared to AD. The 
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Figure 4. Community-wide metrics suggest altered permissivity in niche specificity in colonization in PID individuals. Site codes: (Af) antecubital fossa; 
(N) nares; (Ra) retroauricular crease; (Vf) volar forearm. (A) Mean ± SEM for Shannon diversity is plotted for each patient group (colored as indicated) at all 
sites. (6) Site specificity is measured by the Yue-Clayton theta coefficient as in Figure 3, calculated between moist (antecubital fossa), dry (volar forearm), 
and sebaceous (retroauricular crease) skin sites and nares. (C) Longitudinal stability was assessed by calculating the theta coefficient between two and four 
timepoints for controls (light blue) and STAT3-HIES patients (red). 



decreased S. aureus proportion may relate to the antibiotic pro- 
phylaxis given to reduce S. aurais skin and lung infections, and 
these antibiotics frequently are ineffective against other staphylo- 
coccal species. 

Frequent infections are common in PID patients, and the 
increased frequency or severity of infections is a marker of wors- 
ening disease. Therefore, we examined laboratory studies and 
a skin disease severity score in the STAT3-H1ES patients to identify 
taxa or groups of taxa that correlated with disease severity. We 
calculated Spearman correlation coefficients for the correlation 
between taxonomic relative abundance at the Af, a site of disease 
predilection, to the following: SCORAD, serum IgE, serum lym- 
phocytes, lactate dehydrogenase, and serum eosinophils. Hemo- 
globin was included as a laboratory study that was not expected to 
correlate with disease severity. 

Propionibacterium acnes abundance was significantly inversely 
correlated with skin disease severity (P = 3.7 x 10"°^) (Supple- 
mental Table Sll), while S. aureus was significantly correlated (P = 
0.002). From unsupervised clustering of the correlation coeffi- 
cients obtained for the 34 most abundant taxa (Fig. SB), we iden- 
tified two major clusters. The first, smaller cluster contained taxa in 
addition to 5. aureus, whose relative abundances were positively 
correlated with disease severity. The second, larger cluster con- 
tained potential markers of a healthier skin microbiome in our 
patient population; their relative abundances were inversely cor- 
related with markers of disease severity. In addition to P. acnes, 
this group included the Firmicutes Lactobacillus (P = 0.05) and 
several Proteobacteria, including Stenotrophomonas (P = 0.01), 
Novosphingobium (P = 0.04), and Burkholderia (P = 0.01). In- 
terestingly, while the majority of staphylococci in the STAT3-H1ES 
individuals were not S. aureus, these non-S. aureus staphylococci 
were not correlated with disease severity. 

Fungal characterization of PID 

In addition to bacterial infections, patients with STAT3-H1ES often 
develop fungal infections, including candidiasis and Aspergillus 
pneumonias. Thus, we performed ITSl sequencing in a pilot study 
of STAT3-H1ES patients. At nondisease sites, the Vf and N, we ob- 
served a trend toward increased fungal richness in the STAT3-H1ES 
individuals (Fig. 6A; Supplemental Table SI 2), which would further 
support the hypothesis of a more permissive microbial niche in the 
PID skin. This increase in diversity appeared to involve a reduction 



in the relative abundance of the common skin fungi Malassezia (P < 
0.05, all sites) (Fig. 6B; Supplemental Table S12). Moreover, we 
observed a significant increase in the relative abundance of As- 
pergillus (P < 0.06) and Candida (P < 0.05) at all sites in the STAT3- 
HIES patients (Supplemental Table SI 2; Supplemental Fig. S8). 
Since antifungal medications are commonly administered to 
STAT3-H1ES patients, we examined their medication histories for 
use of antifungals. We observed Aspergillus even in the STAT3-HIES 
patients who had not previously received any antifungal treatment. 
These findings further suggest that the altered immunity that 
modifies the potential bacterial microbiome niche in PID pa- 
tients also applies to fungal skin/mucosal communities and may 
contribute to the clinical observation of increased fungal in- 
fections in this patient population. 

Fungal-bacterial associations in PID 

To evaluate major fungal-bacterial associations, we per- 
formed a strain-by-strain co-occurrence analysis based on 
a Spearman correlation of fungal and bacterial taxonomic rel- 
ative abundances at each site. Generally, patterns of correla- 
tion/anti-correlation of relative abundances differed between 
controls and STAT3-H1ES heatmaps (Supplemental Fig. S9). For 
example, in the Ra, Malassezia globosa was generally correlated 
with most genera of Firmicutes and Actinobacteria in both con- 
trols and STAT3-H1ES patients and Malassezia restricta was anti- 
correlated. In the STAT3-H1ES group, streptococci, coryneforms, 
and S. aureus were uniquely correlated with the presence of Can- 
dida, Aspergillus, Cryptococcus, and others. In the Af and Vf, 
streptococci, staphylococci, and Corynebacterium were signifi- 
cantly correlated with non-Malassezia fungi in STAT3-H1ES in- 
dividuals, while controls exhibited lesser interactions between 
these taxa. 

To assess co-occurrence and community assembly on a pop- 
ulation level, including both fungal and bacterial taxa, we cal- 
culated the C-score for each site, which reflects exclusion patterns 
of different taxa within a niche, i.e., the tendency of each taxa to 
exclude one another from a given niche (Stone and Roberts 1990). 
To impute significance to the score, we compared it against the 
mean score obtained from a null distribution of 10,000 permuta- 
tions drawn from the same control or STAT3-H1ES matrix. Based 
on statistical significance of C-scores, the microbiomes at each 
site in the STAT3-H1ES and control groups were not composed of 
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Figure S. Taxonomic association with clinical features. (A) Pie charts comparing the mean relative 
representation of major staphylococcal species across patient groups at all sites. Site codes: (Af) ante- 
cubital fossa; (N) nares; (Ra) retroauricular crease; (Vf) volar forearm. (6) Correlation of taxonomy w'tth 
clinical markers of disease severity for the antecubital fossa. Unsupervised hierarchical clustering of 
Spearman correlation coefficients for the 34 taxa v\/hose mean abundance was >0.25% across STAT3- 
HIES patients. Correlations are calculated against blood levels of lactate dehydrogenase, lymphocytes, 
IgE, and eosinophils, and a score of skin disease severity (scoring atopic dermatitis [SCORAD]). He- 
moglobin was a control. Reds indicate correlation; blues indicate anti-correlation. 



interacting bacterial and fungal consortia (P > 0.05, STAT3-H1ES Vf 
[P = 0.01]), suggesting no significant pattern to community-wide 
co-occurrence in the skin. While these results suggest on a broader 
scale that communities of the skin do not have obvious assembly 
rules, co-occurrence of specific fungal- 
bacterial partners begin to suggest some 
level of interaction between members in /\ 
the two kingdoms. Further investigation 
is required to determine the extent of 
such an interaction, which may occur at 
a local level that is not reflected by -c 
these broader microbial surveys. | 



microbiota of several groups of PID pa- 
tients compared with classical AD patients 
and healthy controls. Our study offers a 
unique opportunity to study potential al- 
terations in the human microbiome in the 
setting of PID. We focused on rare inher- 
ited disorders that shared the cutaneous 
phenotype of AD-like skin disease, so that 
we could identify differences in the 
microbiota that may reflect the specific 
immunological and genetic features of 
these distinct populations. Our analyses 
demonstrate an increased permissivity in 
the bacterial and fungal colonization of 
the skin of PID patients, dysbiosis in com- 
munity diversity, association of bacterial 
communities with skin disease manifes- 
tations and severity, and unique bacterial- 
fungal taxonomic co-occurrences. 

Based on our observations of de- 
creased site specificity and longitudinal 
microbiome stability, as well as coloni- 
zation by unique taxa, we hypothesize 
that PID increases the host's ability to be 
colonized by atypical microbiota. In- 
terestingly, this increased permissive- 
ness of skin did not result in significant 
increases in microbial diversity or coloni- 
zation by novel phyla. This suggests that 
host-specific constraints on the taxonomic 
diversity of skin microbiota, e.g., site 
physiology, persist despite defective host 
immunity in both the recognition and the 
clearance of cutaneous taxa. While the is- 
sue of bioburden, which may or may not 
be significantly altered in these patient 
populations, is important and remains a challenge, our study fo- 
cused on taxonomic diversity. 

We found that presence of skin disease had a dominant 
signature that appeared to be strongly associated with community 



Discussion 

Pathogens are typically considered as dis- 
crete causative agents in virulence studies, 
yet they exist within a rich milieu of other 
microbial species that can influence 
pathogenicity. Disease could derive from 
dysbiosis of the microbial community 
without an invading pathogen domi- 
nating the community. To begin to in- 
vestigate the complex interrelationships 
between host, microbe, and disease, we 
characterized the bacterial and fungal 
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Figure 6. Characterization of fungal skin communities suggest increased permissivity and elevated 
levels of opportunistic pathogenic fungi. (A) Boxplots showing fungal richness (species observed) of 
STAT3-HIES versus controls at all sites; black bars indicate median. (6) Select ITS1 fungal taxonomic 
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composition. In AD, the presence of skin disease was very strongly 
associated with the presence of S. aureus. PID patients are frequentiy 
treated with antimicrobial therapeutics to manage or prevent 
chronic or recurrent skin and systemic infections, yet the impact of 
systemic antibiotics on intact skin is incompletely understood. We 
determined that while S. aureus was most significantly correlated 
with disease severity in these patients, other taxa such as Clos- 
tridiales and Corynebacterium were correlated as well. In contrast 
to AD patients, the primary staphylococcal species in the PID pa- 
tients were S. epidermidis and S. haemolytkus rather than S. aureus, in 
spite of the presence of disease. In these PID patients, antibiotics, 
including penicillins and trimethoprim-sulfamethoxazole, are ad- 
ministered to reduce S. aureus infections and appear to effectively 
reduce S. aureus relative abundance. Whether the non-S. aureus 
staphylococci significantiy contribute to skin disease is currently 
unknown; in our previous study of atopic dermatitis, we observed 
a significant increase in the relative abundance of S. epidermidis 
concurrent with increases in S. aureus (Kong et al. 2012). A possi- 
bility is that these other staphylococci may share a mutuaUstic or 
commensal relationship with S. aureus by enabling greater re- 
sistance to antimicrobial peptides (Peschel et al. 2001; Sieprawska- 
Lupa et al. 2004; Lai et al. 2007) or antibiotics (Wielders et al. 2001; 
Berglund and Soderquist 2008) or by potentiating the S. aureus toxic 
response as with C. albicans (Peters and Noverr 2013). Alternately, 
these species could take advantage of inflamed skin conditions and 
co-colonize with S. aureus (Nilsson et al. 1998; McCrea et cil. 2000; 
Cho et al. 2001; Williams et al. 2002). 

While PID skin microbiomes shared the reduced diversity 
seen in AD, we observed microbial traits unique to the skin disease 
of these PID patients. We compared PID patients with AD patients 
who also were undergoing tteatment for their recurring skin dis- 
ease. While treatment, particularly with antibiotics, may skew the 
relative representation of different taxa, our results reflect the day- 
to-day state in which these patients remain at risk for infections 
despite prophylactic antibiotic treatment. Of note, our cohort of 
WAS patients had less skin disease and were less likely to be on 
systemic antimicrobials compared with the other PID groups in 
our study. This may have contributed to the observed similarities 
in skin bacterial community diversity (except the Ra) between WAS 
patients and healthy controls. However, the increased bacterial 
community diversity in the Ra of the untreated and treated WAS 
patients compared with healthy controls suggests that host factors, 
i.e., immunodeficiency, are likely important causes of the observed 
skin microbiome differences in our PID patients, even when use of 
antimicrobials is taken into account. 

Interactions with skin commensals are key in developing 
cutaneous immunity (Naik et al. 2012). One possibility is that the 
skin microbes in our patient cohort may alter the activation of 
cutaneous immunity. The markedly elevated IgE levels in these 
patients may be related to allergy or atopy, or may reflect an ab- 
errant response to a colonizing microbe. Our data support the 
hypothesis that the skin in PID patients exhibit different selective 
pressures leading to permissiveness for different organisms such 
that an environmental microbe, 5. marcescens, uniquely colonized 
PID skin. Its longitudinal fluctuations, as well as its ability to col- 
onize both characteristically affected and unaffected sites (Af, Ra, 
and Vf), suggest that S. marcescens may be a ttcinsient taking ad- 
vantage of the ecological permissiveness of the immunodeficient 
skin, as opposed to being pathogenic. Notably, S. marcescens in- 
fections in chronic granulomatous disease patients are typically 
treated with the antibiotics, e.g., trimethoprim-sulfamethoxazole, 
which were used in our PID patient populations. In addition. 



several members of the Clostridiales family also significantly colo- 
nized PID skin, whereas In conttols, their representation was very 
low. 

Since STAT3-HIES patients suffer from oral and pulmonary 
fungal infections and DOCK8 patients experience multiple re- 
calcitrant viral infections, investigating nonbacterial microbiota is 
of Importance. The observation of Increased Candida and Asper- 
gillus species on the skin of STATS -HIES patients may correlate with 
their risk of pulmonary aspergillosis and oral candidiasis with the 
skin serving as a possible reservoir for infection, much as the N is 
a reservoir for S. aureus and Streptococcus pneumoniae infections 
(von Eiff et al. 2001). The presence of Aspergillus did not differ in 
patients on antifungal medications with activity against Aspergillus 
compared with patients without a history of antifungal agents 
with activity against Aspergillus, again highlighting the significant 
skin microbial differences observable in PID patients in the ab- 
sence of use of antimicrobial agents. In addition, coinfection of 
certain C. albicans isolates with S. aureus can result in infectious 
synergism with enhanced toxicity and host Inflammation (Peters 
and Noverr 2013). 

From a genetic and mechanistic perspective, WAS and 
DOCK8 deficiency syndromes likely share many similarities. WAS 
acts downstream from DOCK8 as an effector of CDC42 signaling 
(McGhee and Chatila 2010). WAS patients have cytoskeletal 
dysfunction leading to immune deficiency via defective migra- 
tion or Immune synapse formation by a wide range of immune 
cells, and DOCK8-deflcient patients are also thought to have 
cytoskeletal abnormalities based on known functions of other 
DOCK family members. In contrast, defects in STAT3 signaling, 
a key regulator of many immunologic pathways, impair defensin 
expression and neutrophil recruitment and production, in part 
due to defects in Thl 7 differentiation. A distinction among the 
PID groups was that severe skin disease was not observed in these 
WAS patients. Interestingly, the overall skin microbiomes in 
DOCK8-deficient and STAT3-H1ES patients were more similar to 
each other, perhaps reflecting their phenotypic rather than 
mechanistic convergence. Concurrent host/microbiome tran- 
scriptional profiling will further our understanding of the PID 
immunologic response to shifts in their resident microbiome, and 
functional metabolomlcs will examine the microclimate of PID 
skin and how the microbiota respond to, or in turn, affect the 
host. 

While PIDs are rare monogenic disorders, periods of immune 
suppression are frequent in cancer patients and transplant re- 
cipients. Future studies should investigate whether loss of im- 
munity similarly alters microbial community in these states. 
Many of these patients receive prophylactic antibiotic treatment 
to decrease bacterial load, preventing both infection and micro- 
bial translocation. Further characterization of dysbiosis in im- 
mune-compromised individuals, as well as the antibiotic re- 
sistance reservoir in these chronically treated patients, might 
direct both usage of antibiotic prophylaxis and Infection control. 
Future studies may also explore a role for prebiotics/probiotics to 
maintain microbial diversity and health as preventative ap- 
proaches. 

By Koch's postulates, the relationship between pathogen and 
disease denotes chciracteristics of necessity and sufficiency. How- 
ever, it is increasingly clear that microbial pathogens, such as 
5. aureus, which can asymptomatically colonize human hosts, 
may only be disease-associated in the context of their microbial 
communities and host factors. A corollary to this assertion is that 
dysbiosis of the microbial community can itself function as an 
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etiological agent. General instability and dysbiosis of the skin 
microbiome, as observed in PID patients, may be an underlying 
contributor to their recurrent skin disease. If therapies target not 
only pathogenic agents but also microbial homeostasis, the risks of 
developing antibiotic-resistant bacteria might be reduced. In 
summary, the increased instability and permissiveness of the 
bacterial PID skin microbiome lends insight into host-microbiome 
interactions and its potential contribution to a distinct PlD-asso- 
ciated microbiome. 

Methods 

Subjects 

All samples for this study were obtained according to a protocol 
approved by the Institutional Review Board of NHGRI, NIH. 

Healthy and AD subjects were recruited from the greater Wash- 
ington, DC, metropolitan region. PID patients were recruited from 
patient cohorts evaluated at the Clinical Center, NIH. All subjects 
gave informed consent. Demographic and clinical data were col- 
lected for each subject. All subjects provided medical and medi- 
cation history and underwent dermatologic examination at all 
clinic visits. The diagnoses of PID patients were based on clinical 
history and genetic testing. Assessment of extent of skin disease 
involvement was based on the objective scoring atopic dermatitis 
(SCORAD) (Kunz et al. 1997). Bloodwork included serum IgE (lU/ 
mL; normal range: 0-90), hemoglobin (g/dL; normal range: males 
13.7-17.5, females 11.2-15.7), absolute lymphocyte count (K/(jlL; 
normal range: 1.18-3.74), absolute eosinophil count (K/(jlL; nor- 
mal range: 0.04-0.36), and lactate dehydrogenase (U/L; normal 
range: 1 13-226). A summary of subject characteristics are shown in 
Table 1, with detailed information about sampling, sites analyzed, 
and date sampled in Supplemental Table SI. For Sanger data, data 
for healthy controls and AD patients were obtained from previous 
studies (Grice et al. 2009; Kong et al. 2012). 

Sample collection and processing 

All samples were collected by H.H.K. as described in Grice et al. 
(2009) and Kong et al. (2012). Briefly, skin preparation instructions 
included avoiding bathing and avoiding emollients or antimicro- 
bial soaps or shampoos for 24 h prior to all sampling. Sampling 
sites included the nares (N), retroauricular crease (Ra), antecubital 
fossa (At), and volar forearm (Vf). From a 4-cm^ area, bacterial 
swabs (via Epicentre swabs) and scrapes (via sterile disposable 
surgical blade) were obtained and incubated in enzymatic lysis 
buffer and lysozyme for 30 min at 37°C. For the nares, only swab 
samples were taken. Scrapes and swabs were obtained from the 
other sites as available. Correlation analysis of log proportions 
of genera identified in scrape versus swab sampling confirmed 
the similarity of techniques (R = 0.87, P < 2.2 X 10"^*; data not 
shown). 

Preparation of samples for I6S Sanger sequencing 

Samples were processed as described in Oh et al. (2012). Briefly, 
skin swabs were incubated in enzymatic lysis buffer and lysozyme 
(20 mg/mL) for 30 min at 37°C. Two 5-mm stainless steel beads 
(Qiagen) were added to the solution, placed in a TissueLyser 
(Qiagen), and processed for 2 min at 30 Hz. The standard proto- 
col for the PureLink Genomic DNA kit (Invitrogen) was followed 
for all subsequent steps. Near full-length 16S rRNA genes (16S) 
were amplified from purified genomic DNA using primers 8F 
(5'-AGAGTrTGATCCTGGCTCAG-3') and 1391R (5'-GACGGGCG 
GTGWGTRCA-3 ' ) and were sequenced at the National Institutes of 



Health Intramural Sequencing Center (NISC). Three hundred to 
400 unique sequences were obtained from 16S rRNA ampllcons of 
each sample. 

Preparation of samples for 4S4 sequencing 

Genomic DNA extractions were performed as for the Sanger se- 
quencing. 16S rRNA V1-V3 amplicon libraries were prepared from 
sample DNA using AccuPrime HF Taq (Invitrogen, 12346-086) and 
universal primers flanking variable regions VI (27F, 5'-AGAG 
TTTGATCCTGGCTCAG-3') and V3 (534R, 5'-ATTACCGCGGCT 
GCTGG-3'). For each sample, the universal primers were tagged 
with unique sequences ("barcodes") to allow for multiplexing/ 
demultiplexing (Lennon et al. 2010). For ITSl amplicon library 
preparation, universal primers flanking the ITSl region — adapter + 
18SF (5'-CCTATCCCCTGTGTGCCTTGGCAGTCTCAGGTAAAAG 
TCGTAACAAGGTTTC-3') and 5.8S-1R -i- barcode (5'-GTTCAAAG 
AYTCGATGATTCAC-3')— were used to amplify a select number of 
genomic DNA extractions. PCR products were then purified using 
the Agencourt AMPure XP Kit (A63880) and quantitated using the 
Quant-iT dsDNA high-sensitivity assay kit (Invitrogen, Q33120). 
Approximately equivalent amounts of each PCR product were 
then pooled and purified with a Qiagen MinElute column (Qiagen, 
28004) into 30 (jiL TE prior to sequencing at the NIH Intramural 
Sequencing Center. Amplicon libraries were sequenced on a 454 
GS FLX (Roche) instrument using titanium chemistry. 



16S rRNA amplicon sequence analysis pipeline 

Supplemental Table SI shows samples for which analyzable data 
were obtained, as a certain number of samples failed during the 

sequencing and quality control pipeline. Mothur version 1.26.0 
(Schloss et al. 2009) was used for sequence processing, definition of 
operational taxonomic units (OTUs), and downstream analyses. 
For full-length 16S rRNA amplicon processing, sequence assembly, 
filtering, and alignment were performed as described in Grice et al. 
(2009) and Kong et al. (2012). Briefly, sequences matching the 
human genome were removed (E-value < 0.1) and then aligned 
using the Greengenes NAST (DeSantis et al. 2006) aligner. For VI- 
V3 region 16S rRNA ampllcons, 454 flowgram data were denoised 
using the mothur implementation of PyroNoise (Quince et al. 
2011). The resulting sequences were trimmed of primers and 
barcodes (tolerant of primer errors < 1 and barcode errors <2) and 
filtered for homopolymers £8 and length >200 bp. Sequences 
were aligned in mothur using a SILVA reference database (Pruesse 
et al. 2007). 

Chimera removal was performed with the mothur imple- 
mentation of UCHIME (Edgar et al. 2011) and classified using 
a ribosomal database project naive Bayesian classifier (Wang et al. 
2007). Full-length and V1-V3 16S-rRNA sequences were processed 
in parallel, and to ensure that generally results from either method 
resulted in equivalent data, we calculated the partial Spearman 
correlation coefficient between the relative abundances of the 
dominant taxa (Supplemental Fig. SI), adjusting for multiple pa- 
tient timepoint measurements. 

Following sequence quality control processing, all uncor- 
rected pairwise distances were calculated and operational taxo- 
noiTuc units (OTUs) defined at 97% similarity using average 
neighbor joining. No lane masking was applied. To estimate sam- 
pling saturation, rarefaction curves were generated for each site for 
each sequencing method used (Supplemental Fig. S2A,B). While 
OTUs did not achieve saturation, the rarefaction curves confirmed 
that our sampling gave reasonable coverage to analyze the dom- 
inant members of the bacterial communities. Alpha diversity 
(community evenness and richness: Shannon diversity index. 
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"diversity") and beta diversity (shared community structure/ 
membership: Theta, 8 "similarity" index) were calculated at a 97% 
similarity cutoff in mothur. 

Species-level classification for Staphylococcus and Streptococcus 
was performed using the pplacer algorithm (vl.l.alphal3-0- 
glec7786, (32)) in which each 16S rRNAgene classifying to either 
of these genera are compared with a curated set of reference 
packages created as described in Matsen et al. (2010). Briefly, a 
high-quality 16S rRNA reference database was built from rRNA 
records extracted from RefSeq genomes (as of April 2012) and RDP 
type species sequences (release 10, update 24), and sequences for 
each genus were then extracted from this trusted set and placed 
on a phylogenetlc tree using the taxtastic suite of tools (taxit 
vO.3.1) (Conlan et al. 2012). Sequences classified as Staphylococcus 
or Streptococcus by the RDP naive Bayesian classifier were then 
placed on a phylogenetlc reference tree with "-keep-at-most 1000- 
max-pitches 1000" and taxonomic classifications were generated 
using the guppy program in pplacer, with likelihood cutoffs set 
to 0.60. 

ITS! amplicon sequence analysis pipeline 

Sequences were denoised, trimmed of barcodes and primers, and 
chimeras removed with a pipeline as previously described for skin 
fungal community analyses (Findley et al. 2013). Sequences were 
then classified to the genus level using this database and the 
k-nearest neighbor (knn) implementation in mothur and then 
assigned to phylotypes. Sequences of Malassezia species were 
identified using a custom curated Malassezia reference library as 
described in Findley et al. (2013). 

Statistics 

All data are represented as mean ± SEM unless otherwise indicated. 
Unless otherwise indicated, P-values were adjusted for multiple 
comparisons greater than six using the p. adjust function in R using 
method = "fdr" (Benjamini and Hochberg 1995). Statistical sig- 
nificance was ascribed to an alpha level of the adjusted P-values < 
0.1. Each site was treated as a separate data set based on spatial 
physiological differences between different body niches (Grice 
et al. 2009). All per-sample calculations are available upon request 
as a separate data set. For Ra, statistical comparisons were calcu- 
lated for both age-adjusted and non-age-adjusted values. 

For subjects with multiple measurements (e.g., symmetric 
sites or multiple sampling timepoints), samples were processed 
independently and symmetric sites were combined. Multiple 
sampling timepoints were adjusted for using a mixed effects model 
with a random intercept for multiple timepoints. Processed se- 
quences were then subsampled to 1000 sequences per sample (454 
only), and diversity statistics and proportions were calculated from 
the subsampled data. For alpha and beta diversity statistics, pools 
were subsampled lOOx and resultant values averaged. More than 
five subsamplings yielded a 99.8% correlation with results from 
one subsampling, and five to 100 subsamplings were up to 99.99% 
similar. For longitudinal analyses based on the theta coefficient, we 
calculated the average similarity between two and four timepoints 
within an individual and then averaged the mean longitudinal 
stability of all individuals for either STAT3-H1ES or control pa- 
tients. For spatial (intersite) diversity, we calculated the mean theta 
coefficient between moist (antecubital fossa), dry (volar forearm), 
and sebaceous (retroauricular crease) skin sites and the nares for 
each individual, and then calculated the mean for all individuals 
within a category. For mean theta indices, multiple timepoints 
were averaged prior to calculating a between-category average. Fi- 
nally, the nonparametric Kruskal-Wallis test was used with 6 to 



determine statistically significant differences between microbial 
populations. To identify significcint intercategory comparisons, we 
used a post-hoc multiple comparison test, implemented by the 
kruskalmc test in the pgirmess package in R. For ITSl-derived 
measurements, symmetric site and multiple timepoint measure- 
ments were averaged prior to statistical testing. 

For taxonomy-based analyses, we did not subsample typical 
filtering strategies (e.g., taxa must occur in >25% of individuals 
within a patient category) or remove low-abundance/low-incidence 
taxa that generally result from unequal sampling. Where they 
existed, we used relative abundances obtained from full-length 
Sanger sequencing in place of the relative abundances obtained 
from 454 sequencing. To identify taxa that were differentially 
abundant between controls and PID patients, we employed 
a mixed effects model using category as a factor variable and taxa 
relative abundance as the dependent variable, including a ran- 
dom intercept for multiple patient measurements. P-values were 
generated from the mixed effects model using the glht test (cat- 
egory = "Tukey") in the multcomp package in R. Random Forests 
analysis was implemented by the randomForest package in R on 
the filtered data. C-scores were calculated using a null hypothesis 
of random community assembly calculated against 10,000 ran- 
dom matrices through the vegan and bipartite R packages. For 
correlations with patient clinical metadata, we selected taxa 
whose mean relative abundance across the STAT3-HIES patient 
group was >0.25% and calculated a Spearman correlation co- 
efficient with clinical blood markers IgE, eosinophil, lympho- 
cytes, lactate dehydrogenase, and hemoglobin levels, and 
SCORAD. For heatmaps, we visualized those taxa with mean 
relative abundance >0.25%. 

Data access 

The sequence data from this study have been submitted to NCBl 
BioProject (http://www.ncbi.nlm.nih.gov/bioproject) under ac- 
cession number 46333. Patient and sample metadata have been 
submitted to dbGAp (http://www.ncbi.nlm.nih.gov/gap) under 
accession number phs000266.v3.pl. 
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